#ifndef __common
#define __common
#include"physicalconstants.h"
#include<cmath>
#include "particles.h"
/// 根据给定的 激光 omega，返回等离子体临界密度
inline double ncrit(const double& omega)
{
    return (e_mass * epsilon_0 * omega * omega / e_charge / e_charge);
}

/// 根据给定的等离子体密度，返回等离子体频率
inline double omegape(const double& density)
{
    return sqrt(density * e_charge * e_charge / e_mass /epsilon_0);
}
///根据给定的使用omega归一化的矢势，返回激光强度

inline double get_intensity(const double& scaled_A, const double& omega)
{
    /// wavelength in m
    double lambda = c_light_speed / omega * 2.0 * PI;
    /// 
    double ingauss = (scaled_A * scaled_A) * 1.37 * 1e18 * (um / lambda) * (um / lambda);
    return ingauss * w_cm2;
}
/***
double get_intensity(const double& scaled_A, const double& omega)
{
    return (e_mass * c_light_speed * omega * scaled_A / e_charge * 1e-3);
}
***/
/// 粒子数据结构体，用于MPI封装
struct __Pdata{
    /// 粒子种类
    int sp;
    /// 粒子的id和weight
    long id_wt[2];
    /// 粒子的速度，位置和能量
    double vel_pos_ek[14];
};

/// 粒子种类信息结构体，用于MPI数据传输后的数据恢复
struct __Spinfo{
    /// 粒子种类
    int sp;
    /// 粒子质量
    double mass;
    /// 粒子电荷
    double charge;
};


/// 将粒子信息拷贝到mypdata结构体,用于MPI传输
void copy_to_pdata(__Ptcls& ptcl, __Pdata& mypdata)
{
    mypdata.sp = ptcl.species;
    mypdata.id_wt[0] = ptcl.id;
    mypdata.id_wt[1] = ptcl.weight;
    mypdata.vel_pos_ek[0] = ptcl.velocity.member[0];
    mypdata.vel_pos_ek[1] = ptcl.velocity.member[1];
    mypdata.vel_pos_ek[2] = ptcl.velocity.member[2];
    mypdata.vel_pos_ek[3] = ptcl.old_vel.member[0];
    mypdata.vel_pos_ek[4] = ptcl.old_vel.member[1];
    mypdata.vel_pos_ek[5] = ptcl.old_vel.member[2];
    mypdata.vel_pos_ek[6] = ptcl.position.member[0];
    mypdata.vel_pos_ek[7] = ptcl.position.member[1];
    mypdata.vel_pos_ek[8] = ptcl.old_pos.member[0];
    mypdata.vel_pos_ek[9] = ptcl.old_pos.member[1];
    mypdata.vel_pos_ek[10] = ptcl.energy;
    mypdata.vel_pos_ek[11] = ptcl.lorentz_gamma;
    mypdata.vel_pos_ek[12] = ptcl.mass;
    mypdata.vel_pos_ek[13] = ptcl.charge;
}

/// 拷贝结构体粒子到粒子信息
void copy_to_ptcl(__Pdata& pdata, __Ptcls& ptcl)
{
    
    ptcl.species = pdata.sp;
    ptcl.id = pdata.id_wt[0];
    ptcl.weight = pdata.id_wt[1];
    
    ptcl.velocity.member[0] = pdata.vel_pos_ek[0];
    ptcl.velocity.member[1] = pdata.vel_pos_ek[1];
    ptcl.velocity.member[2] = pdata.vel_pos_ek[2];
    ptcl.old_vel.member[0] = pdata.vel_pos_ek[3];
    ptcl.old_vel.member[1] = pdata.vel_pos_ek[4];
    ptcl.old_vel.member[2] = pdata.vel_pos_ek[5];
    ptcl.position.member[0] = pdata.vel_pos_ek[6];
    ptcl.position.member[1] = pdata.vel_pos_ek[7];
    ptcl.old_pos.member[0] = pdata.vel_pos_ek[8];
    ptcl.old_pos.member[1] = pdata.vel_pos_ek[9];
    ptcl.energy = pdata.vel_pos_ek[10];
    ptcl.lorentz_gamma = pdata.vel_pos_ek[11];
    ptcl.mass = pdata.vel_pos_ek[12];
    ptcl.charge = pdata.vel_pos_ek[13];
}

/// 恢复粒子的种类信息根据spinfo
void set_spinfo(__Ptcls& ptcls, __Spinfo& spinfo)
{
    ptcls.charge = spinfo.charge;
    ptcls.mass = spinfo.mass;
}

/// 定义Pline为一个结构体数组,MPI传输使用
typedef std::vector<__Pdata> __Pline;



inline void shift_data_tags(long data_tags[])
{
    for(int i = 0; i < 16; i ++)
    {
        data_tags[i] += 100;
    }
}


inline double calc_gamma(const __Vect3<double>& vel)
{
    return 1.0 / sqrt(1.0 - (vel * vel) / c_light_speed / c_light_speed);
}


inline void grid_wxy(const double& fracx, const double& fracy, double wx[], double wy[])
{
    //double fx, fy;
    double cfx2 = fracx * fracx;
    double cfy2 = fracy * fracy;
    /***
    if(fracx > 0.5)
    {
        fx = fx - 1.0;
        cfx2 = fx * fx;
        wx[0] = 0.25 + cfx2 + fx;
        wx[1] = 1.5 - 2.0 * cfx2;
        wx[2] = 0.25 + cfx2 - fx;
        wx[3] = 0;
        wx[4] = 0;
    }
    else if(fracx < -0.5)
    {
        fx = 1.0 + fx;
        cfx2 = fx * fx;
        wx[0] = 0;
        wx[1] = 0;
        wx[2] = 0.25 + cfx2 + fx;
        wx[3] = 1.5 - 2.0 * cfx2;
        wx[4] = 0.25 + cfx2 - fx;
    }
    ***/
    //else{
        wx[0] = 0;
        wx[1] = 0.25 + cfx2 + fracx;
        wx[2] = 1.5 - 2.0 * cfx2;
        wx[3] = 0.25 + cfx2 - fracx;
        wx[4] = 0;
    //}

    //if(abs(fracy) <= 0.5)
    //{
        wy[0] = 0;
        wy[1] = 0.25 + cfy2 + fracy;
        wy[2] = 1.5 - 2.0 * cfy2;
        wy[3] = 0.25 + cfy2 - fracy;
        wy[4] = 0;
    /***
    }
    else if(fracy > 0.5)
    {
        fy = fy - 1.0;
        cfy2 = fy * fy;
        wy[0] = 0.25 + cfy2 + fy;
        wy[1] = 1.5 - 2.0 * cfy2;
        wy[2] = 0.25 + cfy2 - fy;
        wy[3] = 0;
        wy[4] = 0;
    }
    else
    {
        fy = fy + 1.0;
        cfy2 = fy * fy;
        wy[0] = 0;
        wy[1] = 0;
        wy[2] = 0.25 + cfy2 + fy;
        wy[3] = 1.5 - 2.0 * cfy2;
        wy[4] = 0.25 + cfy2 - fy;
    }
    ***/
    //cout<<"fracx "<<fracx<<"\t"<<"wx "<<wx[1]<<"\t"<<wx[2]<<"\t"<<wx[3]<<endl;
}



#define _frac2 0.25
#define _xmin -1
#define _ymin -1
#define _xmax 1
#define _ymax 1


#endif
